Modeling, dynamical analysis and numerical simulation of a new 3D cubic Lorenz-like system

Little seems to be considered about the globally exponentially asymptotical stability of parabolic type equilibria and the existence of heteroclinic orbits in the Lorenz-like system with high-order nonlinear terms. To achieve this target, by adding the nonlinear terms yz and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$x^{2}y$$\end{document}x2y to the second equation of the system, this paper introduces the new 3D cubic Lorenz-like system: \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{x}=a(y - x)$$\end{document}x˙=a(y-x), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{y}=b_{1}y+b_{2}yz+b_{3}xz+b_{4}x^{2}y$$\end{document}y˙=b1y+b2yz+b3xz+b4x2y, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\dot{z}= -cz + y^{2}$$\end{document}z˙=-cz+y2, which does not belong to the generalized Lorenz systems family. In addition to giving rise to generic and degenerate pitchfork bifurcation, Hopf bifurcation, hidden Lorenz-like attractors, singularly degenerate heteroclinic cycles with nearby chaotic attractors, etc., one still rigorously proves that not only the parabolic type equilibria \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$S_{x} = \{(x, x, \frac{x^{2}}{c})|x\in \mathbb {R}, c\ne 0\}$$\end{document}Sx={(x,x,x2c)|x∈R,c≠0} are globally exponentially asymptotically stable, but also there exists a pair of symmetrical heteroclinic orbits with respect to the z-axis, as most other Lorenz-like systems. This study may offer new insights into revealing some other novel dynamic characteristics of the Lorenz-like system family.

(2) The existence of a pair of heteroclinic orbits to the origin and a pair of symmetrical equilibria.
Our study outcome not only uncovers the interesting dynamics of the cubic Lorenz-like system family, but also provides a reference on predicting the similar dynamical behaviors of other models, especially the higher dimensional ones.
Therefore, in the ongoing pursuit to determine which experimental conditions may require a more complicated model, the present work may offer characteristics of that 3D cubic Lorenz-like system which may be suitable for comparison with experimental data.
The rest of this paper is arranged as follows. Section "Preliminary" introduces some basic concepts. In Section "The new 3D cubic Lorenz-like system", one formulates a new 3D cubic Lorenz-like system and presents some basic dynamical properties of it, i.e., the Chen-like attractor and Lyapunov exponents, bifurcation analysis, singularly degenerate heteroclinic cycles or normally hyperbolic stable foci with nearby chaotic attractors. Section "Basic behaviors" studies the stability and bifurcation of equilibria by utilizing the center manifold theorem, Routh-Hurwitz criterion, the theory of pitchfork bifurcation, Hopf bifurcation and Lyapunov function. In Section "Existence of heteroclinic orbit", combining concepts of α-limit set, ω-limit set and the theory of Lyapunov function, one proves the existence of heteroclinic orbits. Conclusion remarks are drawn in Section "Conclusions". Preliminary Consider the differential system ẋ = f (x, ξ), where x ∈ R n and ξ ∈ R m are vectors representing phase variables and control parameters respectively. Assume that f is of class C ∞ in R n × R m . Suppose that system has an equilibrium point x = x 0 at ξ = ξ 0 . If at least one eigenvalue of the Jacobian matrix associated with linearized vector field about x 0 is zero or has a zero real part, then x 0 is said to be non-hyperbolic or semi-hyperbolic.
In this paper, system (2) has a line of semi-hyperbolic equilibria S z = {(0, 0, z)|z ∈ R} , given by the z-axis. As the value of z varies, S z are saddles, or foci or nodes normally hyperbolic to the z-axis.
In this paper, we define the set S x = {(x, x, x 2 c )|x ∈ R, c � = 0} to the parabolic type equilibria. Referring to 61 , the generic pitchfork bifurcation is that the restriction of a system to the center manifold is locally topologically equivalent near the bifurcating equilibrium point to one of the following normal forms, ξ = mξ ± ξ 3 . As stated in 18,47,48,51,52,62,63 , for system (2), the degenerate pitchfork bifurcation is defined to be the symmetric bifurcation occurring as the certain parameter crosses the zero value, i.e., c = 0 , due to the line of equilibria existing for c = 0 . The main difference between the generic and degenerate pitchfork bifurcation is that, for c = 0 , the flow of the studied system restricted to the 1D center manifold coincides with the center manifold of the system at the origin, associated with the invariant z-axis, which is filled by equilibrium points if c = 0.
Let the set of points: S (either connected or disconnected) be equilibria of ẋ = f (x, ξ) and D ⊂ R n to be a domain containing S. Let V : D → R be a continuously differentiable function such that V (S) = 0 and then S is globally exponentially asymptotically stable.

The new 3D cubic Lorenz-like system
Based on the Lorenz-like system 1 : one in this section proposes the following 3D autonomous chaotic system: c . Finally, a hidden Lorenz-like attractor coexisting with one saddle in the origin and two stable equilibria is coined based on bifurcation diagrams. The manuscript has been uploaded to the web site: "https://github.com/IjbcThree/Kim", and the interested readers can download it. Set X = (x, y, z) T , system (2) is rewritten as Apparently, system (2) is not topologically equivalent to the generalized Lorenz systems family 2 , and it shows a Chen-like attractor with three Lyapunov exponents: Fig. 1. Furthermore, the chaotic dynamics are examined in the following two cases: In this case, based on Proposition 4.1 and 4.7 in Section "Basic behaviors", the equilibria S ± exist and are stable for b 1 ∈ (0, 1.2210) . This coincides well with the bifurcation diagram in Fig. 2a. In particular, at b 1 = 1 , trajectories of system (2) change from the stable S + to the stable S − , which is a sign to chaos as the ones 14,38 , especially the hidden one illustrated in Fig. 3. Particularly, when 1.180 ≤ b 1 < 1.2210 , there exist chaotic attractors coexisting with stable S ± and the saddle S 0 . While 1.2210 < b 1 < 2.869 , system (2) experiences chaotic behaviors coexisting with unstable S ± and the saddle S 0 . But there are periodic three windows in the chaotic band for 1.576 < b 1 < 1.98 . When 2.869 < b 1 < 4 , there is a period-doubling bifurcation window, which is an important route to chaos and is also similar to its special case [ 1 , Figure 4, p.1888]. ( At this time, from Proposition 4.1, 4.2 and 4.7 in Sectipn "Basic behaviors", the non-isolated or line of semihyperbolic equilibria S z exist for c = 0 , and the equilibria S ± also exist and are stable for c ∈ (3.0234, 10.8865) . For c = 0 , singularly degenerate heteroclinic cycles and normally hyperbolic stable foci S z with nearby chaotic attractors exist, as shown in Fig. 2b, which is in accordance with Figs. 4 and 5, despite a little bit on the parameters b 3 and b 4 . When 0 < c < 0.6 , system (2) undergoes chaotic behaviors coexisting with unstable S ± and the saddle S 0 . While 0.6 < c < 0.2 , there is a period-doubling bifurcation window, foreboding a coming chaos.
(2)  (1) i.e., a special case of system (2), it is a difficult task to detect the hidden Lorenz-like attractors in system (2), which might contribute to the power of nonlinear terms and the number of parameters.
Further, based on the dynamics of S z in Proposition 4.4 in Section "Existence of heteroclinic orbit" and through a detailed numerical study, we may state the following numerical result.
) of each normally hyperbolic saddle S 1 z given in Proposition 4.4 tend to one of the normally hyperbolic stable nodes (resp. foci)  Table 1 when the value of z varies.     www.nature.com/scientificreports/   (1) If z < 0.01701 , the trajectories of system (2) starting from the unstable manifold of S z = (0, 0, z) ultimately toward the normally hyperbolic stable foci S z with z > 0.0175 , forming singularly degenerate heteroclinic cycles, as depicted in Figs. 4a and 8.

Basic behaviors
In this section, the stability and bifurcation of equilibria of system (2) are studied by aid of the center manifold theorem, Routh-Hurwitz criterion, the theory of pitchfork bifurcation and Hopf bifurcation, Lyapunov function and so on.
Firstly, from the algebraic structure of system (2), one easily presents the distribution of equilibrium points in the following proposition.
has three equilibria: S 0 = (0, 0, 0) and a pair of symmetrical equilibria Secondly, for the convenience of determining the stability and bifurcation of equilibria, one has to calculate Jacobian matrix associated vector field of system (2): One can easily calculate the characteristic equations of points of S z , S 0 , S x and S ± : (1) The one of each of S z is (2) The one of S 0 is with 1 = −a , 2 = b 1 and 3 = −c.
(4) The one of S ± is:  Based on the center manifold theorem, one can study the two-parameter family of first-order ordinary differential equations on the center manifold of S 0 : to determine the stability of S 0 near b 1 = 0. Therefore, one arrives at through substituting the expanded expressions of V (u,b 1 ) and S(u,b 1 ): into system (6). Further, the restricted vector field of system (6) on its center manifold is obtained by substituting those expressions in (7) into system (6). Since U(0, 0) = 0, ∂U ∂u u=0,b 1 =0 = 0 and a generic pitchfork bifurcation happens at S 0 according to the pitchfork bifurcation theory 23,64-66 .
(2) When the parameter c crosses the zero value, the family of this vector field crosses this degenerate situation transversally. More precisely speaking, for cb 1 [cb 4 + b 3 + b 2 ] < 0 , the line of equilibria S z existing for c = 0 disappears and equilibria S 0 and S ± appear in system (2).
The proof is over.  www.nature.com/scientificreports/ and a 2D W s loc (S 0 ) containing the z-axis.
Proof The proof is similar to the ones in 37,52,54,56 . One only sketches it. For a > 0 , b 1 > 0 and c > 0 , the eigenvalues of S 0 are 1 = −a < 0 , 2 = b 1 > 0 and 3 = −c < 0 . Thus S 0 has a 2D W s loc containing z-axis, and 1D W u loc (S 0 ) whose appropriate expression is and substituting them into system (2), one obtains the following first-order differential equation In addition, the matrix equation Eq. (10), one has K 1 = 0 , H 2 = 0 and K 2 = (a+b 1 ) 2 2b 1 a 2 . The proposition is thus proved.  Table 2, Table 3 lists the local dynamics of S x , where (2) Moreover, for c = 2a > 0 , b 1 = b 3 = 0 and b 2 = −cb 4 < 0 , each point of S x is globally exponentially asymptotically stable.  www.nature.com/scientificreports/ Proof (1) Firstly, the local stability of points of S z and S x easily follows from the linear analysis and is omitted here.
(2) Secondly, we discuss the global stability of points of S x , i.e., each point of S x is globally exponentially asymptotically stable. For c = 2a > 0 , b 1 = b 3 = 0 and b 2 = −cb 4 < 0 , set the following Lyapunov function: with which yields Namely, points of S x are globally exponentially asymptotically stable. The proof is finished.
(2) If c = 0 , then the dynamics of S 0 are the same to the ones of S z with z = 0 and listed in Table 2.
Proof The local stability of S 0 easily follows from the linear analysis and is omitted here.
In the following Section "Conclusions", one studies the existence of heteroclinic orbits of system (2). For the convenience of discussion in the sequel, the following notations are introduced.
Based on Proposition 5.1, the existence of heteroclinic orbits to S 0 and S ± is derived in the following statement. (13) dV (φ t (q 0 )) dt  (a) Neither homoclinic orbits nor heteroclinic orbits to S + and S − exist in system (2). (b) System (2) has a pair of symmetrical heteroclinic orbits to S 0 and S ± .
Proof a) Firstly, one proves that neither heteroclinic orbits nor homoclinic orbits to S + and S − exist in system (2). Assume by contrary that φ t (q 0 ) is a heteroclinic orbit or a homoclinic orbit to S + and S − , i.e.
Finally, let us show that if system (2) has a second heteroclinic orbit to S 0 and S + , then it coincides with γ + . Suppose φ 1 t (q 0 ) is a solution of system (2) that where s 1 − and s 1 . Therefore, one obtains s 1 − = S 0 and s 1 + = S + , i.e., It follows from Proposition 5.1(ii) that φ 1 t (q 0 ) = γ + . Since system (2) is symmetrical with respect to the z-axis, γ − is another unique heteroclinic orbit to S 0 and S − . Figure 10 verifies the correctness of the theoretical result. Thus proof is completed.

Conclusions
This note reports a new 3D cubic Lorenz-like system, which contains the existing one as special cases and generates rich dynamics, such as generic and degenerate pitchfork bifurcation, Hopf bifurcation, infinitely many singularly degenerate heteroclinic cycles with nearby Chen-like attractors, etc. Using Lyapunov functions, we prove that the parabolic type equilibria are globally exponentially asymptotically stable, and there exists a pair of heteroclinic orbits to the origin and two symmetrical equilibria. In future work, other important dynamics of that system, such as homoclinic orbit, invariant algebraic surface, positively invariant set, the forming mechanism of chaotic attractor and so on, require further analytical descriptions to complete its mathematical treatment. We also hope that the basic ideas and the self-contained approach presented in this paper can be applied to explore other similar chaotic/hyperchaotic systems, i.e., where a = 0 , c, b i , p 1 , p 1 , p 2 , q 1 , q 2 , q 3 ∈ R , i = 1, 2, 3, 4 , etc. Preliminary studies show that there might exist chaotic/hyperchaotic attractors, generic and degenerate pitchfork bifurcation, Hopf bifurcation, singularly degenerate heteroclinic cycles, globally exponentially asymptotically stable parabolic type equilibria, a pair of heteroclinic orbits to the origin and two symmetrical equilibria in system (16)(17)(18). We also guess that the dynamics also exist in higher dimensional analogues. (16)     ẋ = a(y − x), y = b 1 y + b 2 yz + b 3 xz + b 4 x 2 y + p 1 w, z = −cz + y 2 , w = q 1 w + q 2 y,     ẋ = a(y − x), y = b 2 yz + b 3 xz + b 4 x 2 y + p 1 w, z = −cz + y 2 , w = q 1 w + q 3 x,     ẋ = a(y − x), y = b 1 y + b 2 yz + b 3 xz + b 4 x 2 y + p 1 w + p 2 x, z = −cz + y 2 + xy, w = q 1 w + q 2 y + q 3 x,